
source("1.0 MyFunc.R")

Boxplot = function(n, scenario){
    # scenario = c(pi_true,deltad_true,delta_true,op_true)
    
    n.name = paste("n=",n,sep="")
    
    if(all(scenario == c(1,1,1,1))) main = "All correct"     # 1
    if(all(scenario == c(1,1,1,2))) main = "M2, M3 correct"  # 1
    if(all(scenario[1:3] == c(1,1,2))) main = "M2 correct"   # 2
    if(all(scenario[1:3] == c(1,2,1))) main = "M3 correct"   # 2
    if(all(scenario[2:3] == c(2,2))) main = "All wrong"      # 4
    if(scenario[1] == 2)        main = "All wrong" # except for 7-2
    if(all(scenario == c(2,1,1,1))) main = "M1 correct"      # 1

    pdf(paste("./Figures/boxplot_n",n,"_scenario_",
              scenario[1],scenario[4],scenario[2],scenario[3],".pdf",sep=""),
              onefile=FALSE, width = 6.0, height = 6.3)
    
    plot.width   <- getOption("width")
    main = Wrap(main,width = plot.width)
    
    result = Results[,scenario[1],scenario[2],scenario[3],scenario[4],
                     n.name,c("reg","b-ipw","g", "b-tr")]
    plot.data = cbind(c(result),factor(rep(colnames(result),each=nrow(result)),
                                       levels=colnames(result)))
    colnames(plot.data) = c("estimate","estimator")
    
    boxplot(estimate~estimator, data = plot.data, main=main, xaxt='n',
            boxwex=0.5, cex.main=1.5, lwd=1.5,cex.axis=1.5,
            # ylim=range(results[,N,,]))  
            ylim=c(-1,1))  
    axis(1, at=1:4,labels=c("b-reg","b-ipw","g","b-tr"),cex.axis=1.5)
    abline(h=Delta.true,col="red",lwd=1.5)
    
    dev.off()
    
}